clear;

%---------------------------------------------------------
%---------------------------------------------------------
disp([' '; ' ']); 
disp('Assemble records into a single matrix (All)'); 
disp(' ');
count=0;
All.time=-780:0;
q.regions=1:17;
for ctr=1:length(q.regions);
    
    %load data    
    if q.regions(ctr)<=16, 
        R=region_data_all(q.regions(ctr));
        fname=['data/data_',R.name,'.mat'];
        
        %Only use data with magnetics. 
        load(fname,'C');
        T=C.T;
        subset=1:length(T);
    end;
    if q.regions(ctr)==17,       
        %load Juan de Fuca perpendicular gridded analysis
        load data/compile_zones_JF.mat Tall;    
        T=Tall; clear Tall; 
        subset=1:length(T);
        R.name='AT000';
        for ctt=1:length(T); 
            d=distance(T(ctt).lon(1),T(ctt).lat(1),T(ctt).lon,T(ctt).lat)*6.378*10^6*100*2*pi/360;
        end;
    end;
    
    %combine into a single array
    for ct1=1:length(T);
        for ct2=1:2,
            if ct2==1,
                t=T(ct1).t;
                z=T(ct1).iz;
            else,
                t=-fliplr(T(ct1).t);
                z=fliplr(T(ct1).iz);
            end;                            
            if nanmean(z)>0, z=-z; end;
         
            pl=find(t>=All.time(1) & t<=0 & isnan(z)==0 & isnan(t)==0);
            if length(pl)>=0.25*length(All.time),
                pl=find(isnan(z)==0 & isnan(t)==0);
                temp=interp1(t(pl),z(pl),All.time);
                if sum(~isnan(temp))>0.25*length(All.time),
                    count=count+1; 
                    if exist('T(ct1).ridge'),
                        All.loc(count,:)=T(ct1).ridge.rxy;
                    else,
                        All.loc(count,:)=[T(ct1).lon(1) T(ct1).lat(1)];
                    end;
                    All.z(count,:)=temp;
                    All.zm(count,:)=temp-nanmean(temp);
                    
                    %difference and remove outliers
                    nsig=6;
                    temp=diff(temp);
                    temp=temp-nanmean(temp);
                    temp(isnan(temp))=0;
                    temp=temp(1:end-5);  %remove end effects
                    temp=temp-nanmean(temp); 
                    stemp=nanstd(temp);
                    pl=find(temp>nsig*stemp);
                    temp(pl)=nsig*std(temp);
                    pl=find(temp<-nsig*stemp);
                    temp(pl)=-nsig*std(temp);
                    temp=temp-nanmean(temp);
                    
                    %spectral analysis
                    [All.P(count,:) All.freq dummy dummy]=pmtmPH(temp,1,4,0);
                    if count==1, 
                        df=diff(All.freq(1:2)); 
                        pl41=find(All.freq>=1/41-df/2 & All.freq<=1/41+df/2);
                        pl100=find(All.freq>=1/(3*41) & All.freq<=1/(2*41));
                    end;
                    [dummy plm]=max(All.P(count,:));
                    All.peak(count,[1 2])=0;
                    if sum(pl100==plm)==1,
                        All.peak(count,1)=1; 
                    end;                                
                    if sum(pl41==plm)==1,
                        All.peak(count,2)=1; 
                    end;
                    All.region(count)=ctr; 
                    All.S(count)=T(ct1).S(min([length(T(ct1).S) ct2])); 
                    
                    if ctr==17 %correct transect count for structure of JDF data
                        if ct1<7
                            All.transect(count,:)=[q.regions(ctr) subset(ct1) ct2];
                        elseif ct1==7 || ct1==8
                            All.transect(count,:)=[q.regions(ctr),subset(ct1)-3 2];
                        else
                            All.transect(count,:)=[q.regions(ctr),subset(ct1)-2 ct2];
                        end 
                    else %all other ridges
                    All.transect(count,:)=[q.regions(ctr) subset(ct1) ct2];
                    end
                    
                    All.zw(count,:)=temp;
                    All.side(count)=ct2; 
                                        
                    disp([R.name,' ',num2str([q.regions(ctr) ct1 ct2 All.peak(count,:)])]),
                end;
            end;
        end;
    end;
end;

save collect_all.mat All;

